UNCLASSIFIED 


Defense  Technical  Information  Center 
Compilation  Part  Notice 

ADP013639 

TITLE:  Residual  Scale  Dependence  and  PDE  Filters 
DISTRIBUTION:  Approved  for  public  release,  distribution  unlimited 


This  paper  is  part  of  the  following  report: 

TITLE:  DNS/LES  Progress  and  Challenges.  Proceedings  of  the  Third 
AFOSR  International  Conference  on  DNS/LES 

To  order  the  complete  compilation  report,  use:  ADA412801 

The  component  part  is  provided  here  to  allow  users  access  to  individually  authored  sections 
of  proceedings,  annals,  symposia,  etc.  However,  the  component  should  be  considered  within 
the  context  of  the  overall  compilation  report  and  not  as  a stand-alone  technical  report. 

The  following  component  part  numbers  comprise  the  compilation  report: 

ADP013620  thru  ADPO 13707 


UNCLASSIFIED 


RESIDUAL  SCALE  DEPENDENCE  AND  PDE  FILTERS 


G.  PANTELIS 

School  of  Aerospace,  Mechanical  and  Mechatronic  Engineering, 
University  of  Sydney,  NSW,  2006,  Australia. 


1.  Introduction 

In  [1]  and  [2]  the  filter  operation  acting  upon  nonlinear  partial  differential 
equations  (PDE)  is  examined  by  way  of  a PDE  filter  on  a domain  involv- 
ing space/time  and  an  additional  dimension  associated  with  a space  scale 
parameter.  Under  this  approach  it  is  possible  to  obtain  an  estimate  for 
the  error  associated  with  the  equations  satisfied  by  the  filtered  solutions 
of  the  microscopic  scale  PDE  given  any  approximation  of  the  residuals. 
This  provides  a condition  of  consistency.  The  PDE  filter  approach  suggests 
approximations  for  the  residuals  that  are  independent  of  empirical,  or  ar- 
bitrary, parameters. 

Here  the  main  points  of  [1]  and  [2]  are  presented  with  some  additional 
issues  addressed.  An  attempt  is  made  to  remain  within  a setting  of  general 
nonlinear  PDE  systems.  The  filtered  equations  of  reactive  turbulent  flows 
are  presented  as  an  example. 

2.  Macroscopic  Equations 

Let  Q C R3,  T — (0,to)  and  g E I = {0,go),  where  to  > 0 and  0 < go  « 1. 
Define  M C Rn  such  that  M = Q x T x /.  If  SI  C R3  is  bounded 
we  denote  the  boundary  of  H by  dSl.  The  summation  convention  of  re- 
peated upper  and  lower  indices  is  adopted  and  for  convenience  the  fol- 
lowing ranges  for  the  indicated  indices  should  be  assumed  throughout: 
a,  b,c,d  E {1, ...,  n — 2};  i,j,  k,l  E {1, ...,  n — 1};  p,  q,r,s  E {1, ...,  n}. 
Here  n (=  5)  is  the  dimension  of  M and  under  this  convention  ( xa ) = 
(^°)i<a<3  = ( x,y,z ),  xn~l  = t,  xn  = g.  Hence  we  write,  Q has  a coor- 

dinate system  (xa),  QxT  has  a coordinate  system  (xl)  = ( xa , t ) and  M has 
a coordinate  system  ( xp ) = (x' ,g).  The  Greek  indices  are  within  the  range: 
a,  (3,7,6  E IV},  and  are  associated  with  the  dependent  variables,  ua, 
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where  N denotes  the  number  of  dependent  variables.  Throughout  we  shall 
work  with  nondimensional  variables. 

The  scale  parameter  77  is  defined  by 

V = ^/2  (1) 


where  e = l/L  represents  the  ratio  of  the  microscopic  length  scale  l (small- 
est resolvable  characteristic  length  scale)  and  the  macroscopic  length  scale 
of  the  domain  L.  The  parameter  /3  (>  0)  is  chosen  to  control  the  rate  of 
damping  of  fluctuations  in  the  dependent  variables.  In  most  practical  appli- 
cations, e represents  a characteristic  grid  size  for  a (nondimensional)  spatial 
discretization. 

The  set  of  smooth  functions  on  M will  be  denoted  by  F(M).  For  any 
/ G F(M)  the  comma  followed  by  a subscript  ftP  will  indicate  the  partial 
derivative  of  / with  respect  to  the  coordinate  xp.  It  should  be  noted  that 
much  of  the  work  that  follows  could  be  relaxed  by  replacing  the  smoothness 
condition  to  one  of  differentiability  of  some  finite  order.  It  will  be  useful  to 
introduce  a subset  F(M)  C F(M)  which  is  defined  as  follows: 

Definition  1.  Any  / G F(M)  has  the  following  properties: 

(i) /GF(M) 

(ii)  For  any  ua  G F(M),  f = fn (ua , u1*, u*j) . 

Property  (ii)  of  the  definition  states  that  any  member  of  F(M)  can  be 
expressed  in  a functional  form  explicitely  independent  of  xp,  ii“  and 

It  should  be  pointed  out  that  F(M)  is  introduced  here  only  to  avoid  some 
notational  difficulties.  In  [2]  the  ideas  that  follow  are  more  conveniently 
developed  in  the  setting  of  contact  manifolds  where  some  of  the  restrictions 
introduced  here  can  be  avoided. 

For  any  ua  G F{M)  let  tpa  G F(M)  be  defined  by 


where  L is  an  elliptic  differential  operator  on  Q. 
Define  the  vector  field  operator 


u = Hu^de  + (L(»'i)),iaj  + 


with  the  notation 


Q 9 ai  d art  d 

9l>  due*  f~But 


A calculation  shows  that 

8 


(2) 

(3) 

(4) 


/ e F(M) 


(5) 
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where 

W = <pada  + <pffi  + <p%3%  (6) 

The  PDE  filter  is  defined  as  follows: 

Definition  2.  Any  ua  G PDEF{L,fa,ua}  has  the  following  proper- 
ties: 

(i)  ua  G F(M),  ua  G F(n  x T)  and 

Ua\r,= 0 = U01  (7) 

(ii)  L : F(M)  F(M)  is  an  elliptic  operator  on  0 and 


(^-IK  = 0,  (xp)  G M 


(iii)  fa  G F(M)  vanishes  at  r/  = 0,  i.e. 

fVo  = /“(u^uj^ujy)  = 0 


In  Definition  2 the  ua(xa,t)  are  exact  solutions  of  the  equations  that 
describe  the  system  at  the  microscopic  scale.  These  are  expressed  in  the 
form  of  the  PDE  system  (9)  defined  on  fl  x T.  ;,From  the  prescribed 
data  ua\v=o  = ua(xa,t)  we  can  generate  a one  parameter  family  of  fields 
ua(xa,t,ri)  on  M by  integrating  (8)  with  respect  to  the  scale  parameter 
r).  The  operator  L must  be  chosen  such  that  the  integration  of  (8),  with 
respect  to  rj,  will  damp  out  fluctuations  in  each  ua  that  cannot  be  resolved 
on  each  scale  associated  with  the  parameter  rj.  We  will  refer  to  ua(xa,t,r]) 
for  r]  > 0 generated  by  (8)  as  the  filtered  fields  associated  with  solutions 
of  the  system  of  PDE  (9).  The  existence  of  each  ua  G F(fl  x T)  (cou- 
pled with  suitable  boundary  conditions  for  the  filtered  fields  on  90  if  0 
is  bounded)  is  sufficient  to  guarantee  the  existence  of  the  filtered  fields 
ua  G PDEF{LJa,ua}. 

An  effective  filter  in  its  simplest  form  can  be  obtained  using  the  space 


Laplace  operator 


L = 


Scd 


d 2 

dxcdxd 


(10) 


where  6cd  = 1 if  c = d and  5cd  = 0 if  c ^ d.  For  demonstration  purposes 
we  assume  (10)  throughout.  It  should  be  noted  that  here  we  are  employing 
space  filters.  Definition  2 could  be  easily  modified  to  include  space  and/or 
time  filters.  The  simplest  time  only  filter  can  be  obtained  by  setting  L = 
d 2 /dt2. 

If  O = R3  and  each  ua  G F(0  x T)  is  bounded  on  SlxT  then  the 
filtered  fields  obtained  from  the  integration  of  (8)  with  respect  to  rj  subject 
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to  the  initial  data  (7)  yields  the  Gaussian  filter 


ua(xa,t,T))=  f G(xa  — za,r])ua(za,t)dz1dz2dz3 

Jq 


where 


G(xa,r])  = (47t?7)  3/2exp[-6abxaxb/{4ri)\ 

Here  the  Gaussian  filter  is  expressed  in  normalized  form  so  that 


[ G(xa  — za,r])dzldz2dzz  = 
J n 


(11) 

(12) 

(13) 


For  the  case  of  bounded  Q we  must  impose  boundary  conditions  for  the 
filtered  fields  on  dfi.  The  PDE  filter,  as  it  is  defined  here,  will  still  be  ef- 
fective in  damping  out  irresolvable  fluctuations  provided  suitable  boundary 
conditions  can  be  imposed. 

As  t/  is  increased  from  0 each  fa  G F(M)  may  deviate  from  0.  We 
introduce 

fa  = ja+ra  (14) 

for  some  residual  ra  G F(M)  such  that 

r%=o  = 0 (15) 

Let'  ea  G F(M)  be  defined  by 

e°  = (f --L)r°-<,°  (16) 

where  oa  G F(M ) is  given  by 

aa  = (L-  U)fa  (17) 

iFrom  (14),  (16)  and  (17)  is  obtained 

Q 

ea  = (—-L)fa-Wfa  (18) 


If  ua  G PDEF{L,  fa,ua}  then  ipa  — 0 and  the  vector  field  operator  W 
vanishes.  Hence  for  ua  G PDEF{L,  fa,ua} 

(£j  - 0/“  = <=“  (19) 

^From  (9)  and  (15)  we  have 


n,=o  = o 


(20) 
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Consider  first  the  case  14  = R3.  We  see  that  each  /“  satisfies  the  initial 
value  problem  (19)-(20),  where  ?/  replaces  the  traditional  role  taken  by  time 
and  t,  appears  only  as  a parameter.  If  each  ea  E F(M)  is  bounded  on  M 
the  solution  to  this  initial  value  problem  can  be  obtained  in  the  explicit 
form  (see  Lecture  8 [3]) 

fa(xa,t,ri)  = r [ G(xa  - za,rj  -Z)eQ(za,t,Odz1dz2dz3d(,  (21) 

Jo  J n 

Thus  for  consistency  we  need  to  generate  residuals  ra  such  that,  through 
the  identity  (16),  each  e°  is  rendered  sufficiently  small. 

If  0 has  a boundary,  Oil.  then  we  can  assume  that  given  some  pre- 
scription of  ua  and  rQ  (and/or  their  spatial  gradients)  on  the  boundary 
dO,  we  can  set  fa\gn  = 0.  This  along  with  the  system  (19)-(20)  defines  an 
initial  boundary  value  problem  satisfied  by  each  /“.  For  an  estimate  of  the 
consistency  error  we  can  replace  the  expression  (21)  by  some  inequality. 
For  instance  in  terms  of  the  I/2(f4)-norm,  ||  • ||,  we  can  obtain,  for  each 
a E {1, ...,  N},  T]  6 I and  t E T, 

ii  n-,x)  ii<  /Vhon-^-e  (22> 

Jo 

for  some  A > 0.  The  inequality  also  holds  for  the  case  A = 0. 

If  ea  = 0 on  M then  the  residuals  ra  are  known  exactly.  Thus,  on  M, 
the  filtered  fields  of  the  solutions  of  the  PDE  (9)  and  the  exact  residuals 
rQ  E F(M)  satisfy  the  system 

/“  + ra  = 0 (23) 

A 

(vr L)ra  — <7°  = 0 (24) 

or] 

where  oa  is  given  by  (17).  These  are  the  exact  macroscopic  equations  for  the 
filtered  fields  ua  E PDEF {L,  fa,  uQ}.  Written  in  this  form  it  is  clear  where 
the  difficulties  arise  in  constructing  subgrid  scale  models  for  general  nonlin- 
ear PDE  systems  (9).  In  application  one  desires  the  solution  computed  only 
on  some  slice  M|7?=const.  The  presence  of  the  term  dra  jdr\  makes  the  solu- 
tion of  (24)  impractical  and  some  approximation  for  the  residuals  needs  to 
be  introduced.  The  task  is  to  generate  ra  such  that  eQ  is  minimzed  in  (16). 
The  relationship  (21)  (or  (22))  gives  an  estimate  of  how  close  (23)  is  satis- 
fied by  the  filtered  fields  ua  E PDEF{L,  /“ ,ua}  given  any  approximation 
of  the  residuals  ra. 

3.  Approximation  of  the  Residuals 

We  consider  an  approximation  to  (24),  presented  in  [1]  and  [2],  that  is  valid 
for  any  slice  M\v=const  for  0 < 77  < 770 • Assume  ra  E F(M),  such  that 
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ra|^=o  = 0,  satisfy  on  M 

(-  - L)ra  -aa  = 0 (25) 

V 

for  some  ua  E F(M),  not  necessarily  in  PDEF{L)  fa,  £tQ},  and  cr"  given 
by  (17).  A calculation  based  on  (15),  (16)  and  (25)  gives 


e = 


dra 

dr) 


1 n d2va 

5K-r«U)-lV 


l,=£,  {€(0,<»)  (26) 


If  d2ra/dr]2  can  be  suitably  bounded  then  we  have  ea  = 0(r/).  ^From  (21) 
(or  (22))  it  follows  that  fa  — 0(rj2)  for  ua  E PDEF{L,  fa,ua}. 

In  general  any  ua  E PDEF{L,  /“,  ua}  will  not  satisfy  (23)  exactly  on 
M if  the  residuals  are  not  exact,  i.e.  they  do  not  satisfy  (24).  In  application 
one  generates  ua  E F(M ) by  enforcing  (23)  and  introducing  some  approx- 
imation for  the  residuals  ra  (for  example  obtained  from  (25)).  In  such  a 
case  each  ua  can  only  approximate  members  of  PDEF{L,  fa,  ua  \. 


4.  Approximation  of  the  Filtered  Fields 

Let  ua  E F(M)  satisfy  (i)  and  (iii)  of  Definition  2 and  be  generated  from 
(23)  given  some  approximation  of  the  residuals  ra  (for  example  through 
(25)).  Let  ua  E PDEF{L,  fa,  uQ}.  Since  ua  is  not  a member  of  PDEF{L , 
ua}  there  exists  <pa  E F(M)  satisfying  (2)  and  we  have 

(|j -£)(««  -««)=*>“  (27) 

Consider  first  the  case  J7  is  unbounded  (i.e.  f2  = R3).  If  ipa  E F(M ) is 
bounded  on  M then  we  obtain 

(ua  — ua)(xa,t,rj)  = [ [ G(xa  — za,r]  — £)ipa(za,t,£)dz1dz2dz3d£  (28) 
Jo  J n 

For  the  case  where  fi  has  a boundary  dO.  we  assume  that  ua\en  = ua|an-  As 
before,  in  terms  of  the  L2(0)-norm  we  can  obtain,  for  each  a E {1, ...,  N}, 
rj  E I and  t E T, 


||  (ua  -ua)(;V)  ll<  r II  <paU)  II  (29) 

Jo 

for  some  A > 0.  The  inequality  also  holds  for  the  case  A = 0.  Both  (28)  and 
(29)  suggest  that  convergence  of  the  approximations  of  the  filtered  fields 
rests  upon  the  boundedness  of  (pa. 
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A calculation  based  on  (16),  (17)  and  (23)  gives 

Wfa  = ~ea  (30) 

If  the  residuals  ra  ^ire  exact,  then  ea  = 0 and  Wfa  = 0.  For  this  to 
hold  for  arbitrary  fa  we  must  have  W = 0 (the  zero  operator).  Hence 
ipa  = 0 and  ua  € PDEF{L,  f'\  ua}.  If  each  ra  is  not  exact  then  each  ua, 
generated  from  (23),  will  not  be  a filtered  field  and  each  tpa  will  not  vanish 
everywhere  on  M.  Suppose  that  approximations  for  the  residuals,  ra , can 
be  found  rendering  each  ea  = 0(S(p)),  where  8 is  an  order  function  with 
respect  to  r]  such  that  limn^8(r])  = 0.  The  problem  is  then  transformed 
into  one  of  establishing  the  rate  at  which  the  operator  W tends  to  the  zero 
operator  as  rj  — > 0. 

5.  Example:  Reactive  Flows 

The  ideas  of  the  previous  sections  are  developed  in  the  context  of  general 
nonlinear  PDE.  For  demonstration  purposes  we  assume  an  incompressible 
fluid.  The  fa  are  associated  with  the  three  fluid  momentum  equations  and 
/n_1  is  associated  with  the  continuity  or  fluid  mass  conservation  equation. 
The  system  is  augmented  with  fn~l+A  (yl  = 1,...,Q)  associated  with  Q 
mass  balance  equations  for  the  Q chemical  components.  In  this  case  we 
have  N = n — 1 + Q and 

ua  = va , un~l=p , un-1+A  = ujA,  (A  = 1,...,Q)  (31) 

where  va  (1  < a < 3)  correspond  to  the  fluid  velocity  components,  p the 
fluid  pressure  and  u)A  (A  = 1, ...,  Q)  are  the  mass  fractions  of  the  chemical 
components.  We  can  write  the  functions  fa  € F(M)  corresponding  to  the 
equations  of  motion  of  an  incompressible  fluid  and  chemical  component 
mass  balance  as 

fa  = v“  + (vbva  + 8abp  - 8bcvyR.e)fi 
fll~l  - vbjb  (32) 

fn-l+A  = UA  + (vbuA  _ KSbcu>A):b  + ZA 

where  Re  is  the  Reynolds  number,  k is  a coefficient  associated  with  molec- 
ular diffusion  of  the  chemical  components  in  the  fluid  medium  and  ^A  = 
are  the  chemical  source  terms.  It  should  be  mentioned  that 
to  regard  (9)  under  the  prescription  (32)  as  an  exact  description  of  the 
system  at  p = 0 is  not  entirely  true  since  the  fluid  viscous  and  chemical 
diffusion  terms  are  only  approximate  models  of  nonlinear  effects  ocurring 
at  even  smaller  scales.  For  the  sake  of  this  demonstration  we  shall  regard 
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this  formulation  as  representing  the  exact  microscopic  equations  at  77  = 0. 
We  should  note  that  the  viscous  and  diffusion  terms  are  linear  and  do  not 
contribute  to  the  residuals  under  filtering. 

A calculation  based  on  (17)  and  (32)  leads  to 

<7“  = 2 

a""1  = 0 (33) 

<y"~1+A  = 2S‘d(v^Ad),i  + SA 


where 

= (L  - U)£a  (34) 

is  left  in  the  generic  form  since  we  have  not  specified  the  source  terms 
as  explicit  functions  of  the  chemical  components.  We  can  assume  that  the 
source  terms  £A  G F(M). 

The  identities  (33)  provide  the  source  terms  that  appear  in  the  system 
(23)- (24)  (or  the  approximate  system  (23),  (25)).  Since  an~l  = 0 we  can  set 
rn_1  = 0 in  (16)  and  hence  the  continuity  equation  will  be  invariant  under 
filtering,  i.e.  /n_1  = 0 on  M.  It  is  seen  that  each  cr°,  for  the  incompressible 
case,  depend  only  on  the  velocities  and  their  spatial  partial  derivatives. 
The  corresponding  residuals,  ra,  will  model  the  influence  of  the  residual 
stress/strain  of  the  turbulent  fluid. 

The  source  terms  an~1+A  can  be  decomposed  into  two  parts 


Tn-1-M 


n-l+A  , n-l+A 
udisp  ' ureac 


(35) 


where  ^isp+A  associated  with  dispersion  and  cr™eal+A  is  associated  with 
the  reaction  kinetics.  In  the  absence  of  reaction  source  terms  (i.e.  when 
£a  — 0)  each  cr”ea).+'4  vanish  and  the  residuals,  rn~l+A , will  model  the 
influence  of  the  dispersion  of  the  chemical  components  in  the  flow  field. 
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